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Symbiotes can be transmitted from parents to offspring or horizontally from unrelated 
hosts or the environment. A key question is whether symbiote transmission is similar 
enough to Mendelian gene transmission to generate and maintain coevolutionary 
associations between host and symbiote genes. Recent papers come to opposite 
conclusions, with some suggesting that any horizontal transmission eliminates genetic 
association. These studies are hard to compare owing to arbitrary differences in modeling 
approach, parameter values, and assumptions about selection. I show that associations 
between host and symbiote genes (extra-genomic associations) can be described by the 
same dynamic model as conventional linkage disequilibria between genes in the same 
genome. Thus, covariance between host and symbiote genomes depends on population 
history, geographic structure, selection, and co-transmission rate, just as covariance 
between genes within a genome. The conclusion that horizontal transmission rapidly 
erodes extra-genomic associations is equivalent to the conclusion that recombination 
rapidly erodes associations between genes within a genome. The conclusion is correct 
in the absence of population structure or selection. However, population structure can 
maintain spatial associations between host and symbiote traits, and non-additive selection 
(interspecific epistasis) can generate covariances between host and symbiote genotypes. 
These results can also be applied to cultural or other non-genetic traits. This work 
contributes to a growing consensus that genomic, symbiotic, and gene-culture evolution 
can be analyzed under a common theoretical framework. In terms of coevolutionary 
potential, symbiotes can be viewed as lying on a continuum between the intimacy of 
genes and the indifference of casually co-occurring species. 



Keywords: vertical transmission, horizontal transmission, symbiosis, gene-culture, coevolution, dual-inheritance 
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1. INTRODUCTION 

The view that organism phenotypes can be described in terms of 
a dichotomy between inherited genes and non-inherited environ- 
mental factors has been an enormously useful simplification in 
the development of quantitative genetics and evolutionary theory 
(Lynch and Walsh, 1998; Futuyma, 2009). However, several fields 
now present important opportunities to understand the preva- 
lence and importance of additional influences. For example, stud- 
ies of gene-culture coevolution (Feldman and Zhivotovsky, 1992; 
Henrich et al., 2008) and symbiosis (Bright and Bulgheresi, 2010; 
Gilbert et al., 2010) have emphasized the roles of factors with a 
mixture of horizontal and vertical transmission on development 
and evolution. Formation of intimate symbioses has contributed 
to major transitions in eukaryotic evolution and community ecol- 
ogy (Maynard Smith and Szathmary, 1995; Selosse et al., 2006), 
while gene-culture coevolution can promote speciation in learn- 
ing animals (Vallin and Qvarnstrom, 201 1) and is associated with 
the emergence of Homo sapiens as a global ecosystem engineer 
(Maynard Smith and Szathmary, 1995; Vitousek et al., 1997). 

The degree to which symbiotic associations and cultural traits 
are passed vertically from parent to offspring vs. horizontally 
between unrelated individuals can be critical in determining the 



evolutionary trajectories of traits affecting cooperation, resource 
use, and the functional integration of systems (symbiotic or 
social) more inclusive than the individual organism. Vertical 
transmission is thought to promote partner fidelity feedbacks 
between genes (Frank, 1994). That is, when host and symbiote 
genes maintain a statistical association across generations, the 
evolution of stable mutualism is likely (Frank, 1994; Doebeli and 
Knowlton, 1998; Fletcher and Doebeli, 2009; Wyatt et al, 2013). 
When vertical transmission is perfect, the evolutionary dynamics 
converge on the dynamics of coevolution between genes within 
the same genome. 

Several studies have sought to extend the tools and concepts 
of population genetics to analyze host-symbiote and gene-culture 
coevolution. For example efforts have been very successful in 
advancing understanding and manipulation of incompatibility- 
inducing endosymbiotes (Turelli, 1994; Hoffmann et al., 2011). 
Analyses of gene-culture coevolution have been illuminating but 
controversial when applied to human culture (Cavalli-Sforza 
and Feldman, 1981; Boyd and Richerson, 1985; Richerson and 
Boyd, 2005; Ackland et al, 2007; Claidiere and Andre, 2012; 
Houkes, 2012). Even genetic models of speciation can be extended 
to include interactions between host and symbiote genomes 
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(Brucker and Bordenstein, 2012, 2013). Important technical 
advances include the recognition that gene-culture covariance 
and interspecific disequilibrium (covariance between host and 
symbiote genes) resemble conventional linkage disequilibrium 
and cytonuclear disequilibrium (Feldman and Cavalli-Sforza, 
1984; Sanchez et al., 2000). Further, the notion of "general- 
ized epistasis" (Feldman and Cavalli-Sforza, 1984; Feldman and 
Zhivotovsky, 1992) emphasizes the potential for non-additive 
interactions between genotypes and cultural traits to determine 
the direction of evolution. 

These results have inspired a generalized concept of "non- 
genetic inheritance" to include symbiotes, cultural traits, and 
potentially other non-genetic or epigenetic traits with some possi- 
bility of both vertical and horizontal transmission (Bonduriansky 
and Day, 2009). Proponents of this expanded view of inheritance 
have generally concluded that there are important similarities 
between models of non-genetic inheritance and classic genetic 
models, and that extra-genomic traits can be significant factors in 
evolution (Feldman and Zhivotovsky, 1992; Jablonka and Lamb, 
1998; Bonduriansky and Day, 2009; Day and Bonduriansky, 
2011; Odling-Smee et al., 2013). However, attempts to extend 
classical population and quantitative genetics theory to accom- 
modate non-genetic inheritance have been relatively complex, 
making it difficult to establish general principles (Feldman 
and Zhivotovsky, 1992; Santure and Spencer, 2006; Tal et al., 
2010; Johannes and Colome-Tatche, 2011; Bonduriansky, 2012). 
Moreover, the general importance of non-genetic traits in evo- 
lution remains debated (Jablonka and Lamb, 1998; Haig, 2007; 
Dickins and Dickins, 2008; Bonduriansky, 2012). 

The opposing conclusions from recent theoretical studies most 
likely arise from differing premises, assumptions, and param- 
eter spaces. Early work on gene-culture coevolution adapted 
classical population genetic models to assumptions about cul- 
tural transmission, and tended to emphasize effects of selection 
on cultural traits (Feldman and Cavalli-Sforza, 1984; Boyd and 
Richerson, 1985; Feldman and Cavalli-Sforza, 1989; Feldman 
and Zhivotovsky, 1992). Studies of genomic imprinting have 
generally used a quantitative genetic framework emphasizing 
phenotype distributions in families (Santure and Spencer, 2006; 
Tal et al, 2010; Johannes and Colome-Tatche, 2011). Day and 
Bonduriansky (2011) used the Price equation to develop expres- 
sions for phenotypic change owing to selection, genetic inheri- 
tance, and non-genetic inheritance. Their approach appears to 
be quite general, but as a trade-off for generality, it does not 
immediately reveal the importance of any particular mechanism 
of non-genetic inheritance, nor any particular measure of coevo- 
lutionary association. In contrast, Brandvain et al. (2011) made 
an explicit population genetic model focused on how mater- 
nal transmission of a symbiote affects interspecific covariance 
(disequilibrium) in a panmictic population with no selection. 

Under those conditions, Brandvain et al. (2011) showed that 
genetic covariances (interspecific disequUibria) between neutral 
organelle and symbiote alleles decay rapidly with even a little 
horizontal transmission in a panmictic population. They con- 
cluded that imperfect vertical transmission would leave negligible 
statistical signature in molecular marker data and that the poten- 
tial for interactions between host and symbiote genes to affect 



evolution would be limited. Qualitatively, the same can be said for 
conventional nuclear genes or cytoplasmic and nuclear genomes; 
linkage disequUibria (covariances) decay rapidly in panmictic 
populations (Lewontin, 1974; Hartl and Clark, 1997). However, 
conditions favoring persistent covariance between nuclear genes 
or between cytoplasmic and nuclear genes are common in nature 
(Arnold, 1993; Hewitt, 2001; Zapata et al., 2001; Laurie et al, 
2007). Moreover epistasis is widely acknowledged as an impor- 
tant factor in genome evolution and speciation (Wolf et al., 
2000; Coyne and Orr, 2004; Petkov et al., 2005; Weinreich et al, 
2005; Muir and Moyle, 2009; Rohlfs et al, 2010), and Drown 
et al. (2013) showed conditions under which interspecific epis- 
tasis would favor an evolutionary transition from horizontal to 
vertical transmission in a host-symbiote system. 

The apparently unresolved question is whether associations 
between host and symbiote genes have fundamentally different 
dynamics from associations between genes in the same genome 
(as implied by Brandvain et al, 2011), or whether fundamen- 
tal dynamic similarities trump a superficial distinction between 
genetic and non-genetic inheritance (Day and Bonduriansky, 
2011). Here, I attempt to address this question by deriving a 
few very simple relationships, and then illustrating some of their 
implications with simulations of a few biologically interesting 
scenarios. My conclusion might not be surprising to theoretical 
population geneticists; extra-genomic traits like symbiotes can be 
studied by extending the tools of classical population genetics 
(Cavalli-Sforza and Feldman, 1981; Boyd and Richerson, 1985; 
Sanchez et al., 2000; Day and Bonduriansky, 2011). However, 
my analysis might be helpful for microbial ecologists. I shed 
some new light on the early work on associations (disequUib- 
ria) between classical alleles and non-genetic traits (Feldman 
and Cavalli-Sforza, 1984) and contradict the inferred conclusions 
(though not the mathematical results) of Brandvain et al. (2011). 
My results help illustrate the value of a unified view of inheritance 
and coevolution (Day and Bonduriansky, 2011), with symbiosis 
falling on a continuum of intimacy including genes, organelles, 
free-living species, and even environmental factors (Lewontin, 
1983; Keeling and Archibald, 2008; Odling-Smee et al, 2013). 

2. METHODS 

Brandvain et al. (2011) made two claims against the evolution- 
ary relevance of associations between host organelle and symbiote 
alleles (interspecific disequilibrium). First, high rates of vertical 
transmission of symbiotes do not result in high levels of genetic 
association. Second (consequently), host-symbiote interactions 
will not respond to natural selection on non-additive fitness 
effects (interspecific epistasis). The basis of these claims is a model 
showing rapid decay of covariances between organelle and sym- 
biote genomes in a single randomly mating population with no 
selection. Therefore, I first evaluate the veracity of their result 
using a simpler, more traditional population genetics model. 

The model follows earlier work generalizing the traditional 
measure of association between genes in the same genome (link- 
age disequilibrium) to extra-genomic associations (Feldman and 
CavaUi-Sforza, 1984; Feldman and Zhivotovsky, 1992), but I pro- 
vide a new result, a closed-form solution for the evolutionary 
dynamics of extra-genomic associations. To put a finer point on 
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the question of whether extra-genomic associations can be infor- 
mative about the transmission process, I explicitly compare the 
dynamics of symbiote-organelle covariance with the dynamics of 
symbiote-nuclear covariance. To my knowledge, this has not been 
done previously. 

Second, I evaluate the effects of population structure on pat- 
terns of extra-genomic covariance. Population structure is well 
known to promote within-genome covariance (linkage disequi- 
librium), but the impact of subdivision and gene flow on the 
dynamics of extra-genomic traits has rarely been investigated 
(Ackland et al., 2007). I derive new analytical results for equi- 
librium covariance in a simple admixture model (Asmussen and 
Arnold, 1991) and then show numerical results for a more com- 
plicated hybrid zone model modified from Dakin (2006). 

Finally, I test whether epistatic interactions between host and 
symbiote genes can affect how a population responds to selec- 
tion. I illustrate the earlier results of Feldman and Cavalli-Sforza 
(1989), which demonstrate an evolutionary response to selection 
on epistatic effects between a gene and a cultural trait. I show con- 
ditions under which host-symbiote genetic covariance is built up 
by selection despite the tendency for horizontal transmission to 
reduce covariance. I also evaluate effects of hitchhiking between 
host and symbiote genes in a single population and across a 
hybrid zone. 

3. RESULTS 

3.1. COVARIANCE BETWEEN SYMBIOTE AND ORGANELLE GENOMES 

To derive the expected covariance between organelle (cytoplas- 
mic) and symbiote genomes, consider a system where each 
genome has two types (haplotypes or haplotype groups): C and 
c for the cytoplasmic organelle, and S and s for the symbiote. If Xy 
is the frequency of individuals with cytotype and symbiote j, and 
pj and pj are the marginal frequencies of cytotype i and symbiote 
j, the covariance (cyto- symbiote disequilibrium, D) is 



D = xcs- peps- 



(1) 



Equation ( 1 ) is the covariance between two binary random vari- 
ables, and can be applied to phenotypes, genes, symbiotic states 
or any pair of binary categorical variables (Clark, 1984; Feldman 
and Cavalli-Sforza, 1984). Statistical association or disequilib- 
rium can also be expressed as a correlation (r = D/^pspspcpc)- 
If the probability of vertical (maternal) transmission of the sym- 
biote is V, the set of frequencies x = (xcs, xcs, x^s, x^) is expected 
to change across generations according to the transition matrix in 
Table 1. Therefore, the expected cyto-symbiote covariance in the 
next generation is 



D' = x'cs - peps 

= xcsv + (xcs + xcs) (1 -v)ps- peps 

Substitute pc = xcs + xa 



(2) 



Table 1 | Transition matrix for joint symbiote and organelle 
genotypes. 
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And the general solution is Df = Dgv'. This novel result is directly 
comparable to the classic result for covariance between nuclear 
genes Dt = Dqc' , where c is the probability of cosegregation 
(1— recombination rate) (Lewontin, 1974). 

If the symbiosis is not obligate and exclusive, i.e., if a host can 
have zero symbiotes or more than one symbiote genotype, the 
result is similar. Let xcs be the frequency of individuals with cyto- 
type C and symbiote genotype S, regardless of what other sym- 
biote genotypes are present. Then x'^^^ = xc5v(l — hs) + pchs, 
where hs is the probability of gaining a symbiote with genotype 
S by horizontal transfer (from the environment or another host). 
Let D = Xcs ~ pcpSt where ps is the fraction of hosts with sym- 
biote S (rather than the fraction of symbiotes with genotype S). 
In this case 



Dr = Do[v(l-hs)Y . 



(3) 



That is, when the ecological relationship between host and 
symbiote is facultative, the effective recombination rate is 
1 — v(l — /jj), which is greater than the obligate case (1 — v). 

3.2. COVARIANCE BETWEEN SYMBIOTE AND NUCLEAR GENES 

For the covariance between a symbiote genotype and an allele 
at a nuclear locus (with alleles A and a), the nuclear-symbiote 
covariance is D = xaas + ^-X^AaS ~ pApS (Asmussen et al., 1987). 
The relevant host-symbiote genotype frequencies in the next 
generation are given by 

^AAS = ^ (^XaaS + 2*AaS^ pA + - v) ps (^Xaa + 2^Aa^ PA 



and 



'^AaS 



^aaSpA^ 



+ (1 - V) ps ( XAApa + -XAa + XaapA 



The nuclear-symbiote covariance in the next generation D' 



■AAS 



AaS ■ 



■ pAps simplifies to 



D' = xcsv- pcPsv 
= Dv 



D' 



1 

D-v 
2 



(4) 
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And the general solution is Dt = Dq (^v)*. A pertinent spe- 
cial case is perfect maternal transmission (v = 1), which elimi- 
nates the distinction between symbiote and organelle. Covariance 
between a maternally transmitted organelle and a nuclear gene 
(cyto-nuclear disequilibrium) decays by one half each generation 
{Dt = Do (I) ), just like covariance between unlinked nuclear 
genes (Asmussen et al, 1987). 

Even though all covariances decay rapidly in a panmictic 
population (Lewontin, 1974; Sanchez et al., 2000; Brandvain 
et al., 201 1), cyto-symbiote covariance (Equation 2) decays much 
more slowly (by a factor of 2) than nuclear- symbiote covariance 
(Equation 4) whenever there is non-zero maternal transmission 
(Figure 1). The key difference between host-symbiote transmis- 
sion and conventional gene or organelle transmission is that 
completely free recombination corresponds to a cosegregation 
probability of j , whereas completely horizontal transmission cor- 
responds to a vertical transmission probability of zero. Thus, 
although conventional linkage disequilibrium between unlinked 
loci is expected to decay no faster than by | each generation, host 
symbiote covariance can drop to zero in one generation if there 
is no tendency for parent-offspring transmission. Aside from this 
quantitative difference, the population genetic principles describ- 
ing the relationship between nuclear and cytoplasmic genomes 
(Asmussen et al., 1987; Arnold et al, 1988; Arnold, 1993) can be 
extended to the relationship between hosts and symbiotes. This 
insight underscores a fundamental continuity between nuclear 
genes, organelles, and symbiotes. 

3.3. HOST-SYMBIOTE COVARIANCE IN STRUCTURED POPULATIONS 

3.3. 1. Continent-island structure 

Admixture or immigration between populations with divergent 
allele frequencies is a common cause of within-genome associa- 
tion (conventional linkage disequilibrium) in nature (Conner and 
Hartl, 2004). To assess the effect of this kind of population struc- 
ture on host-symbiote genetic associations, first consider a simple 
continent-island admixture model (Asmussen and Arnold, 1991) 







0.20 






D 

.00 0.10 






o 


0 5 10 15 20 
Generation 


FIGURE 1 1 Decay of CYto-symbiote covariance (solid line. Equation 2) 
and nuclear-symbiote covariance (dashed line. Equation 4) for vertical 
transmission probability v = 0.90. 



in which a focal population receives a proportion m of immi- 
grants from two or more divergent source populations (Figure 2). 
Let the frequency of immigrants with cytoplasmic allele C and 
symbiote S be itcs- After immigration, the covariance between 
cytoplasmic and symbiote genotypes (cyto-symbiote covariance 
Dcs) is 



Dcs = (1 - m) xcs + mncs - peps 

= v{l-m) Dcs + mDm + cov(pc, Ps) 



(5) 



where D,„ is the covariance among immigrants. The averages {pc 
and ps) and covariance of the allele frequencies are over immi- 
grant and resident born sets with weights m and 1 — m. Exactly 
the same dynamics can be derived for conventional nuclear or 
cyto-nuclear linkage disequilibria (Asmussen and Arnold, 1991). 
If we assume the source populations are themselves unchanged by 
gene flow, so that D^ and the allele frequencies are constant, this 
difference equation has a non-zero equilibrium at 



Dcs 



mD„ 



1 



v(l — m) 

Likewise for nuclear-symbiote covariance. 



Dns 



mD„ 



1 - jv{l - m) 



(6) 



(7) 



Setting V = 1 in Equation (7) recovers Asmussen and Arnold's 
(1991) solution for cyto-nuclear disequilibrium. Thus, inter- 
specific genetic covariances can be maintained by immigration. 
Moreover, Dcs > Dms', the expected cyto-symbiote covariance is 
greater than the expected nuclear-symbiote covariance at equilib- 
rium between immigration and recombination. 

3.3.2. Stepping-stone structure 

To investigate the effects of more complicated population struc- 
ture on cyto-symbiote and nuclear-symbiote covariances, I added 
a symbiote to the stepping stone model analyzed thoroughly for 
cyto-nuclear covariance by Dakin (2006). The basic framework 
is a line of demes between two source populations. Dispersal 
occurs only between adjacent demes and mating is random within 
demes (Kimura and Weiss, 1964; Goodisman et al., 1998; Dakin, 
2006). This framework can be modified to include some other 
classic models as special cases, including the continent-island 
(above) and two-population intermixture models (Nei and Li, 
1973; Asmussen and Arnold, 1991). 

First, I determined how to incorporate a symbiote with verti- 
cal transmission probability v into Dakin's (2006) deterministic 
recursion equations (Appendix Equation 9). I then iterated the 
system of equations to illustrate transient and equilibrium pat- 
terns. I followed Dakin (2006) by initializing ten stepping stones, 
the first five (nearest source 1) fixed for A, C, and S, and the 
second five (nearest source 2) fixed for a, c, and s. This simu- 
lates secondary contact. Each time step began with symmetrical 
dispersal between neighboring demes, followed by random mat- 
ing within demes. In one set of simulations, I followed Dakin 
(2006) in assuming infinite source populations unaffected by gene 
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FIGURE 2 I Intermixture model used to illustrate how covariance 
between a cytoplasmic genotype Cand symbiote genotype Scan be 
maintained by gene flow (Equations 5-7). Relative contribution of source 
population 1 is represented by f. xcsa) is the frequency of individuals with 



cytoplasmic genotype C and symbiote genotype S in Source /, and ncs is the 
frequency of that genotypic combination in the immigrant pool, x^j is the 
frequency of that genotypic combination in native-born residents in the 
intermixed population, and m is the immigration rate. 



flow (remaining fixed for their respective alleles for all time). 
In a second set, I allowed gene flow into the source popula- 
tions, assuming they received half the number of immigrants 
(because they have only one rather than two neighbors). I report 
the covariances among loci in the adults (after dispersal, before 
mating); covariances in the offspring are lower in magnitude (by 
factors of V, ^, and j) but the differences between cyto-symbiote, 
nuclear-symbiote, and cyto-nuclear associations are greater. 

Dynamics of the stepping stone model over the first 100 gener- 
ations of contact are illustrated in Figure 3. Over this time frame, 
the difference between the infinite source and finite source models 
was negligible. With moderate immigration (m = 0.10 shown) 
the central demes (at the initial contact front) show transient 
covariances many times higher than their equilibrium values, as 
shown previously for cyto-nuclear covariances (Asmussen and 
Arnold, 1991). Most important, with high but imperfect verti- 
cal transmission (v = 0.90 shown), the cyto-symbiote association 
is greater and shows larger changes in time and space than the 
nuclear-symbiote or cyto-nuclear associations. 

Allele and genotype frequencies were very close to their 
asymptotic values within a few hundred generations. For the 
finite source model, all covariances approached zero, as expected 
(Asmussen and Arnold, 1991). Asymptotic results for the infi- 
nite source model (after 1000 generations) presumed to reflect 
equilibrium are summarized in Figure 4. Cyto-nuclear results for 
the allelic and genotypic covariances agree with those of Dakin 
(2006). Cyto-symbiote associations are many times higher, as 
expected from equations 2 and 4 for a high probability of vertical 
transmission (v = 0.90 shown). 

To check the effects of genetic drift on these predictions, I 
also modeled the system with finite populations. At each genera- 
tion in each deme, I sampled N = 1000 adults from the genotype 



distribution in focal deme with probability I — m and from each 
neighboring deme with probability m/2. Then I sampled N = 
1000 offspring given the adult genotype frequencies and trans- 
mission probabilities. Results for this stochastic version were 
similar to the deterministic results for asymptotic covariances 
(Figure 4). Drift tends to increase genetic associations, but this 
is somewhat masked by its effects on allele frequencies. The cor- 
relation (r = D/^pspspcpc) is more sensitive than the covariance 
(Figure 4) because r accounts for local allele frequencies. 

3.4. GENERALIZED EPISTASIS AND SELECTION: CAHLE AS 
SYMBIOTES 

Feldman and Cavalli-Sforza (1984) pointed out that the covari- 
ance between genes and a culturally transmitted trait could be 
described by a disequilibrium (D or r) in exactly the same way as 
the standard covariance between genes, and Sanchez et al. (2000) 
observed that gene-culture covariance and host-symbiote covari- 
ance were formally analogous to cytonuclear disequilibrium (all 
being covariances between binary random variables). My analysis 
shows that this similarity is deeper than the summary statistics; 
the evolutionary associations in host-symbiote and gene-culture 
systems follow the same kind of dynamics as Mendelian genetic 
systems. It follows that evolutionary relationships between host 
and symbiote genes should respond to epistatic selection (con- 
tra Brandvain et al, 201 1). That is, the notion of epistasis can be 
generalized to include interactions between traits with different 
inheritance mechanisms (Feldman and Zhivotovsky, 1992), just 
as the notion of coevolution can be extended to include interac- 
tions between genes in the same genome (Lovell and Robertson, 
2010). 

To underscore this fundamental correspondence, consider the 
coevolutionary model of human lactose absorption and dairy 
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FIGURE 3 I Dynamics of covariances over the first 100 generations of secondary contact with immigration rate m = 0.10 and vertical transmission 
probability v = 0.90. Only the middle six demes are shown; demes 1 , 2, 9, and 10 have essentially flat lines at this scale. 



farming investigated by Feldman and Cavalli-Sforza (1989). Their 
model included an autosomal (diploid) gene affecting lactase 
activity, and milk use as a binary cultural trait. This is a model 
of interspecific "generalized epistasis" because the fitness effect of 
a genotype depends on the presence or absence of the interaction 
between humans and cattle (Feldman and Cavalli-Sforza, 1989; 
Feldman and Zhivotovsky, 1992). The change in frequency of 
the absorption allele depended strongly on the degree of vertical 
transmission of the cultural trait (the probability that children of 
milk users would become milk users). The model applies equally 
well if we characterize dairy farming as a symbiotic association 



between humans and cattle. If cattle or farms tended to be passed 
from parent to offspring, the gene-culture covariance describ- 
ing an association between human genotype and behavior (milk 
use) is exactly the same as host-symbiote covariance describing 
an association between human genotype and symbiote (cattle) 
presence. 

I implemented the Feldman and Cavalli-Sforza (FCS) model 
as follows (see the online appendix and Feldman and Cavalli- 
Sforza, 1989 for mathematical and coding details). The first step 
is random mating with imperfect mother-offspring transmission 
of the symbiote ( [7 vs. u for present or absent) and Mendelian 
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transmission of an autosomal diploid locus determining whether 
a host benefits or suffers a net fitness cost from the symbiosis 
(Table!). In the second step, following Feldman and Cavalli- 
Sforza (1989), offspring have a chance of picking up the symbiote 
via horizontal transmission. That is, each individual without the 
symbiote has probability fpu of becoming infected (pu is the 
frequency of infection/symbiosis among other hosts in the pop- 
ulation). The third step is selection according to Table 2, and 
the resulting genotypes mate at random to produce the next 
generation. 




CM 

O 
O 



oo 
o 
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o 
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o 
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FIGURE 4 I Asymptotic covariances after 1000 generations of secondary 
contact with immigration rate m = 0.10 and vertical transmission 
probability v= 0.90. Circles show averages of 100 replicates of the 
stochastic version of the model (1000 individuals per deme) and lines show 
the deterministic (no drift) expectations. Cyto-symbiote associations (black 
circles, solid line) are several times greater than nuclear-symbiote (open 
circles, dashed line) and cyto-nuclear (gray) associations. The correlation (r^) 
is affected more by drift. 



Table 2 | Fitness of host genotypes participating in the symbiosis (U) 
or not (u) in the Feldman and Cavalli-Sforza (FCS) model. 

Host genoype 



AA 



Aa 



Infected (U) 
Uninfected (u) 



1 + Si 
1 



1 + Si 
1 



1 -S2 

1 



If both S; and S2 > 0, there is sign epistasis because the effect of the EGT is 
negative or positive depending on the host genotype. 



To illustrate that selection can operate on interspecific epistatic 
effects despite imperfect vertical transmission, I present numeri- 
cal results for a range of vertical transmission rates given selection 
si = 0.05 and S2 = 0.15, and horizontal transmission rate / = 
0.5. These values are among those used by Feldman and Cavalli- 
Sforza (1989) in their discussion of the coevolution of lactase 
and dairying, and I chose them by trial and error for the sim- 
ple purpose of providing a counterexample to the conjecture that 
imperfect transmission wiU make interspecific epistatic effects 
unresponsive to selection (Brandvain et al., 2011). 

Perfect vertical transmission is not required. Rather, epistatic 
fitness effects can drive coevolution between host and symbiote 
genes if the vertical transmission rate is above some threshold 
value determined by the strength and mode of selection and the 
probability of horizontal transmission (Figures). In the specific 
example depicted (Figures), there is a strong deleterious effect 
of the symbiosis (or cultural trait) on the "wild-type" genotype 
(52 = 0.2) relative to its benefit for the mutant (si = 0.05), that 
is, the additive effect of the symbiote is negative when the A host 
allele is rare. Nonetheless, the symbiosis can spread (along with 
the A allele) if vertical transmission is strong enough to main- 
tain an association between the symbiote and the host allele that 
makes it advantageous (v > 0.76 in this case, with/ = 0.5). 

Feldman and Cavalli-Sforza (1989) allowed transmission (v 
and /) to vary with genotype. The resulting dynamics differ 
only quantitatively from the simplified version used here. Day 
and Bonduriansky (2011) analyzed a version of the FCS model 
including self-learning. Individuals could (with some probability) 
choose the cultural practice that best suited their genotype (anal- 
ogous to partner choice in symbiosis). They solved for a threshold 
strength of selection above which milk use increases when rare. 
This threshold was lowered by self-learning. High rates of vertical 
transmission lowered the threshold when self-learning was weak, 
but increased it when self-learning was strong. 

3.5. HITCHHIKING EFFECTS 

To evaluate the potential for selection on extra-genomic traits to 
drive changes in organelle genotype frequencies via hitchhiking, I 
evaluated the correlated response of organelle genotype frequen- 
cies to selection in the basic model, the FCS model, and the hybrid 
zone model. The effect of an advantageous, maternally trans- 
mitted symbiote on a neutral cytoplasmic marker can be found 
for the basic model (Equation 2 and Table 1) by adding a fit- 
ness advantage of 1 -F s for individuals with symbiote S. Then the 
frequency of the cytoplasmic allele C changes as 



Pc 



Xcsil+S)+Xcs 

W 



pc + 



sD 

W 



(8) 



where W = I + sps is mean fitness. 

To evaluate hitchhiking effects of a non-obligate symbiote, 
I added a maternally inherited cytoplasmic marker to the FCS 
model. I numerically iterated the FCS model to illustrate hitch- 
hiking effects for three kinds of conditions. 

First, for comparison to the standard case of a neutral marker 
linked to an advantageous mutation (Maynard Smith and Haigh, 
1974; Barton, 2000), I used a single panmictic population with 
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a uniformly advantageous symbiote in the FCS framework (si = 
—$2 = s). The symbiote was initially in complete disequilibrium 
with a cytoplasmic marker. I.e., there were two kinds of individu- 
als: CU and cu. Each CU produced CU offspring with probability 
V and Cm offspring with probability I — v. Cu and cu offspring 



were converted to CU and cU with probability fpu, where pjj 
was the frequency of individuals with the symbiote at that time. 
Individuals with the symbiote had relative fitness 1 + s. For 
comparison to the standard Maynard-Smith and Haigh (MSH) 
model, I calculated the expected final frequency of the cytoplas- 
mic marker (pc.oo) in terms of initial frequencies of C and 17 and 
selection (s) from Barton (2000) Equation (1): pc.oo = pc,o + 
(l~Pc,o)i'yo' approximate effective recombination 

rate r = 1 — v(l — /) from Equation (3). This MSH expecta- 
tion is not perfectly identical to the cyto-symbiote hitchhiking 
model (there are differences in ploidy, initial conditions, and the 
"effective recombination rate" is only an approximation), but the 
magnitude of hitchhiking should be comparable if the claim is 
true that the FCS model produces evolutionary dynamics for 
host-symbiote gene pairs that are comparable to dynamics for 
conventional genes. 

As expected based on standard population genetics and 
Equations (2, 3, and 8), the hitchhiking effect of a spreading sym- 
biote depends on the strength of selection, the co-transmission 
rate, and the initial conditions (Table 3). The examples in 
Figure 6 illustrate the similarity between cyto-symbiote hitch- 
hiking (or extra-genomic hitchhiking in general) and the hitch- 
hiking effect of one nuclear locus on another (Maynard Smith 
and Haigh, 1974; Barton, 2000). Graphs for other parameter 
combinations are available as Appendix A. 



Table 3 | Cytoplasmic marker frequencies after 100 generations of 
cyto-symbiote hitchhiking in the FCS model for various vertical 
transmission rates (v) and fitness effects. 



V s=0.10 


0.15 


0.20 


Po = 0.001 






0.90 0.006 


0.038 


0.095 


0.95 0.039 


0.100 


0.164 


0.99 0.082 


0.151 


0.215 


Po = 0.50 (admixture) 






0.90 0.711 


0.777 


0.818 


0.95 0.752 


0.805 


0.838 


0.99 0.775 


0.821 


0.849 


Contact zone midpoint 






0.90 0.627 


0.680 


0.712 


0.95 0.671 


0.709 


0.734 


0.99 0.692 


0.725 


0.747 


Contact zone with epistasis 






0.90 0.591 


0.673 


0.725 


0.95 0.655 


0.705 


0.735 


0.99 0.668 


0.708 


0.736 



The first case (po = 0.001) represents an initially rare symbiote in a panmictic 
population, and the second (po = 0.50) represents 1:1 admixture (the cytoplas- 
mic marker was initially perfectly associated with the advantageous symbiote). 
The third and fourth cases represent contact zones (mean allele frequencies of 
demes 50 and 51 out of 100 total stepping stones). In the first three cases, the 
symbiote is universally advantageous fc = Sj = —S2). In the final, epistatic case 

S2= jS,. 
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FIGURE 5 I Coevolution of a "host" gene and presence or absence of a 
cultural/symbiotic trait in the Feldman and Cavalli-Sforza (FCS) model 
for a single population. Panel (A) shows prevalence of hosts (e.g., 
humans) with the symbiote (e.g., dairy cattle) over time for six different 
vertical transmission rates (see legend in B). Blue and black lines approach 
equilibrium values between 0.5 and 1.0 (Note logarithmic y-axis). (B) shows 
concurrent changes in frequency of the allele affecting fitness of hosts with 
and without the symbiotic relationship (e.g., a lactose absorption allele). (C) 
Shows the correlation (interspecific disequilibrium) between the host allele 
and symbiote presence. 
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FIGURE 6 I Cyto-symbiote hitchhiking effect over 100 generations 
with a universally advantageous symbiote in a panmictic 
population. Panel (A) slnows dynamics of symbiote prevalence and 
cytoplasmic marker frequency for an initially rare symbiote (bottom pair 
of lines) and for 1:1 admixture (upper pair) under the basic (CS) model 
with selection s=0.15 and vertical transmission i/= 0.9025. (B) Shows 
final cytoplasmic marker frequencies (after 100 generations) as a 
function of the initial frequency (of both cytoplasmic marker and 



symbiote) for the basic CS model in comparison to the Maynard Smith 
and Haig (MSH) approximation for a pair of nuclear loci with s = 0.15 
and r = 1 — 1/ = 0.0975. (CD) Show analogous results from the 
Feldman and Cavalli-Sforza (PCS) model with a universally advantageous 
symbiote (si=0.15, S2 = — 0.15), vertical transmission probability 
v=0.95 and horizontal transmission probability f=0.05. Here I used a 
recombination rate 1 — i/(l — f) = 0.0975 for the MSH model, making 
the two sets of results as comparable as possible. 



Second, to evaluate the dynamics of secondary contact with 
a universally advantageous symbiote, I used the stepping-stone 
model described above with inheritance and selection according 
to the FCS model after dispersal. I used 100 stepping stones and 
iterated the deterministic difference equations for 100 generations 
after secondary contact for several combinations of transmission 
rates (v and/) and selection strengths (si = —$2 = s). 

Finally, I repeated the simulations of secondary contact 
with selection according to the FCS model with sign epistasis 
(Weinreich et al., 2005). That is the effect of the symbiote is pos- 
itive or negative according the host genotype, as in the original 
formulation (Feldman and Cavalli-Sforza, 1989). For the exam- 
ples depicted, I varied the strength of selection but kept the ratio 
constant (52 = Si/2). 

Cyto-symbiote hitchhiking can cause significant differential 
introgression of cytoplasmic DNA (Figure?, Table 3). The spa- 
tial cline for the allele frequency of a neutral nuclear marker is 
unaffected by the spread of a symbiote across a hybrid zone, 
but there is a significant spatial displacement of the cline for 
a neutral cytoplasmic marker owing to the slow decay of cyto- 
symbiote covariance. As it does for intragenomic hitchhiking 



(Barton, 2000), spatial structure reduces the overall hitchhik- 
ing effect (Table 3) by slowing the response to selection and 
allowing more time for recombination. Introgression is generally 
slower with epistasis (Table 3), in part because the overall selective 
advantage of the symbiote is reduced when S2 > — Si (Table 3), 
and in part because the advantage of the symbiote depends on co- 
introgression of the epistatic nuclear allele, producing a steeper, 
more slowly advancing wave front (Figure 7). 

4. DISCUSSION 

The simple models presented here extend previous theoretical 
work on nuclear, cyto-nuclear, and extra-genomic covariances 
(Nei and Li, 1973; Asmussen and Arnold, 1991; Feldman and 
Zhivotovsky, 1992; Goodisman et al., 1998; Sanchez et al., 2000; 
Dakin, 2006; Brandvain et al., 2011). My results illustrate a very 
important conceptual and mathematical correspondence between 
the dynamics of associations between genes in different genomes 
and associations between genes in the same genome. Although 
there have been doubts about the stability of coevolutionary asso- 
ciations between host and symbiote genes with less than 100% 
vertical transmission (Brandvain et al, 2011), such concerns are 
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FIGURE 7 1 Cyto-symbiote hitchhiking effect after 100 generations of 
the Feldman and Cavalli-Sforza (FCS) model in a secondary contact 
zone with 100 demes in a stepping stone pattern (only the central 20 
demes are shown). Selection favoring the symbiote was Si =0.15, vertical 
transmission probability v = 0.95 and horizontal transmission probability 
t = 0.05. Panel (A) shows results for a universally advantageous symbiote 
with no epistasis (S2 = — si). (B) Shows results for sign epistasis with 

S2 = jSi. 



no more or less valid for extra-genomic associations than for asso- 
ciations between genes within a genome. This inference can be 
extracted from the general frameworks of dual- inheritance theory 
(Cavalli-Sforza and Feldman, 1981; Feldman and Cavalli-Sforza, 
1984; Boyd and Richerson, 1985; Feldman and Zhivotovsky, 
1992) or non-genetic inheritance (Bonduriansky and Day, 2009; 
Day and Bonduriansky, 2011). My analysis uses more specific 
models to put a finer point on the idea that coevolution between 
genomes and within genomes share fundamental similarities. 

The study of covariance between genes in symbiotic species is 
important because assortment between mutually beneficial geno- 
types is a key to the evolutionary origin and stability of mutualism 
(Frank, 1994; Fletcher and Doebeh, 2009; Wyatt et al, 2013). 
Vertical transmission does not guarantee the evolution of mutual- 
ism; persistent conflicts are sometimes evident within genomes or 
between host and organelle or endosymbiote genomes (Burt and 
Trivers, 2006). However, cross-generation covariance between 
host and symbiote genes can strongly affect the direction of 
evolution. 

My analysis shows that extra-genomic covariances follow the 
same dynamics as conventional linkage disequilibria. In situations 
where conventional linkage disequilibria tend to persist (e.g., 
structured populations), extra-genomic covariances can also be 
expected if there is some level of vertical transmission. Therefore, 



extra-genomic covariances might often be informative about rates 
of vertical vs. horizontal transmission. In addition, the potential 
for non-additive fitness effects to generate significant evolution- 
ary change at the level of symbiotic or cultural systems can be 
similar to the potential for epistasis to affect change at the level of 
conventional genomes. 

4.1. SIMILARITIES AND DIFFERENCES BETWEEN SYMBIOTE 
TRANSMISSION, CULTURAL TRANSMISSION, AND GENE 
TRANSMISSION 

The derivations of Equations (2-7) expose a striking dynamic 
similarity of within-genome and between-genome covariances. 
Disequilibria (covariances) are expected to decay as exponen- 
tial functions of time. The consequence in a large, panmictic 
population is that we expect covariances to become negligible 
after a few generations (Lewontin, 1974; Sanchez et al., 2000; 
Brandvain et al, 2011). However, even in such an idealized 
population, the rate of decay depends on the probability of co- 
transmission. Therefore, maternally transmitted symbiotes will 
tend to have higher covariance with maternally transmitted genes 
(e.g., cytoplasmic genomes) than with nuclear genes (Figure 1). 

Two differences between intra-genomic and extra-genomic 
associations can modify the exponential decay pattern. First, zero 
linkage corresponds to a cotransmission probability of 0.5 for 
a pair genes within a genome, whereas zero vertical transmis- 
sion corresponds to a cotransmission probability of 0.0 for host 
and symbiote genes. This means the decay of extra-genomic 
covariance can be instantaneous (when v = 0). Second, if the 
symbiotic relationship is not obligate, some offspring can be with- 
out the symbiote (i.e., they neither inherited it from a parent 
nor picked it up from someone else or the environment). This 
effectively increases the decay of extra-genomic covariance and 
causes a systematic tendency for the symbiosis to be lost in the 
absence of selection or an environmental source (Feldman and 
Cavalh-Sforza, 1989; Lipsitch et al, 1995; Sanchez et al, 2000). 

4.2. USING EXTRA-GENOMIC COVARIANCE TO MAKE INFERENCES 
ABOUT VERTICAL TRANSMISSION 

Given the rapid decay of covariances in panmictic populations, 
even for physically linked genes, lack of evident covariances 
cannot be interpreted as evidence against vertical transmis- 
sion (Lewontin, 1974; Sanchez et al, 2000; Brandvain et al., 
2011). However, when present, covariances can be informa- 
tive. My analysis shows how even imperfect maternal transmis- 
sion of symbiotes generates greater cyto-symbiote association 
than cytonuclear or nuclear-symbiote association (Figure 1). This 
result suggests that extra-genomic covariance can be used to make 
inferences about rates of vertical transmission in structured pop- 
ulations. Cytonuclear covariance and covariance between nuclear 
genes with known linkage relationships can be used as bench- 
marks for the detectability of vertical transmission of symbiotes. 
In particular, maternally inherited symbiotes are expected to have 
greater covariance with maternally inherited genes than with 
biparentally inherited genes (by a factor of 2' in the simple case 
of neutral decay in a panmictic population). Maternal trans- 
mission of bacteria might be common in animals (Funkliouser 
and Bordenstein, 2013). Nuclear and cytonuclear covariances are 
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frequently observed in nature, particularly in zones of admixture 
or hybridization (Arnold, 1993; Hewitt, 2001; Zapata et al, 2001; 
Laurie et al., 2007). Interspecific covariances should be common 
in those same situations for symbiotes with substantial vertical 
transmission. 

Recently, Thierry et al. (2011) documented extra-genomic 
covariance between whitefly genotypes and prevalence of several 
endosymbiotic bacteria in a zone of recent contact between native 
and introduced whiteflies. These endosymbiotes are maternally 
transmitted and showed strong correspondence between mtDNA 
and symbiotypes in hybrids (Thierry et al, 2011). However, the 
symbiotes have varying degrees of horizontal transmission and 
only one of three is obligate (Ahmed et al., 2013). Despite imper- 
fect vertical transmission, host mtDNA and symbiote genotypes 
covary in the admixture zone (Thierry et al., 2011) and at a 
regional scale (Ahmed et al., 2013). These observations from 
nature indicate that it will be feasible to test for patterns of 
extra-genomic covariance consistent with maternal transmission 
in many wild systems. In the whitefly system, it will be interesting 
to compare patterns of nuclear gene flow with those of mtDNA, 
primary endosymbiotes, and secondary endosymbiotes. 

4.3. IMPLICATIONS FOR SYMBIOSIS AND EXTRA-GENOMIC 
COEVOLUTION 

Selection on non-genetic traits such as symbioses or learned 
behaviors can produce hitchhiking effects and extra-genomic 
covariances even when vertical transmission is less than perfect. 
Even non-additive interactions (generalized epistasis) between 
genes in different genomes can affect the direction of evolution. 
This has been shown before, particularly in models of gene- 
culture evolution (Feldman and Cavalli-Sforza, 1989; Ackland 
et al, 2007; Day and Bonduriansky, 2011), but the continu- 
ity of coevolution within genomes, between genomes, and even 
between genes and culture is underappreciated (Brandvain et al, 
2011). 

All forms of coevolution depend on consistent associations 
between partners to produce reciprocal responses to selection, 
but intimate symbioses are often thought of as special. Recent 
enthusiasm over the evolutionary ecology of human-microbe 
interactions has raised broad questions about whether plants 
and animals should be conceptualized as "holobionts" or "meta- 
organisms", including hosts and associated microbial communi- 
ties as units or levels of organization in evolution and ecology 
(Rosenberg et al, 2007; Doolittle and Zhaxybayeva, 2010). The 
answer depends, in part, on whether coevolutionary dynamics 
between genes in different genomes (host and microbe) are sim- 
ilar to dynamics between genes in the same genome. The decay 
of associations between genomes in a panmictic population with 
imperfect vertical transmission might be interpreted as leaving lit- 
tle opportunity for non-additive interactions between genomes 
to drive evolution (Brandvain et al., 201 1). However, my analysis 
and several previous studies illustrate how associations between 
genomes are affected by the same factors and follow the same 
kinds of dynamics as associations between genes within a genome. 
Population structure can generate and maintain associations that 
might affect local coevolution (Doebeli and Knowlton, 1998; 
Thompson, 2005; Wyatt et al, 2013) and generalized epistasis can 



generate associations and determine the outcome of coevolution 
(Feldman and Cavalli-Sforza, 1984; Feldman and Zhivotovsky, 
1992; Drown et al., 2013). Interspecific epistasis can even con- 
tribute to host speciation (Brucker and Bordenstein, 2012), for 
example by impacting mating behavior (Sharon et al., 2010) or 
hybrid viability (Brucker and Bordenstein, 2013). It does not nec- 
essarily follow that host-symbiote systems should be conceptual- 
ized as "meta-organisms," but the theoretical continuity between 
coevolution of genes within genomes and between genomes is 
encouraging for further synthesis between evolutionary genetics 
and evolutionary ecology. 
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APPENDIX 

APPENDIX A: EXPANDED DESCRIPTION OF METHODS FOR SIMULATING EXTRA-GENOMIC TRAIT TRANSMISSION 
Host-symbiont cava ria nee in structured populations 

To investigate the effects of population structure on cyto-symbiont and nuclear-symbiont disequilibria, I added a symbiont to the step- 
ping stone model analyzed for cyto-nuclear disequilibria by Dakin (2006). Within each deme, random mating and inheritance of the 
nuclear locus A, maternally transmitted cytoplasmic locus C and imperfectly transmitted symbiont S follows the set of deterministic 
recursions given in Equation (9). 

The expected frequencies after dispersal are given by the weighted averages of adjacent demes such that a fraction tn of each deme 
is composed of immigrants coming in equal proportions from the two adjacent demes (Kimura and Weiss, 1964; Goodisman et al, 
1998; Dakin, 2006). I followed Dakin (2006) by initializing ten stepping stones, the first five (nearest source 1) fixed for A, C, and S, 
and the second five (nearest source 2) fixed for a, c, and s. This simulates secondary contact. Each time step began with symmetrical 
dispersal between neighboring demes, followed by random mating within demes. In one set of simulations, I followed Dakin (2006) 
in assuming infinite source populations unaffected by gene flow (remaining fixed for their respective alleles for all time). In a second 
set, I allowed gene flow into the source populations, assuming they received half the number of immigrants (because they have only 
one rather than two neighbors). I report the disequilibria among loci in the adults (after dispersal, before mating); disequilibria in 
the offspring are lower in magnitude (by factors of v, j, and j) but the differences between cyto-symbiont, nuclear-symbiont, and 
cyto-nuclear associations are greater. 

I iterated the system of equations to illustrate transient and equilibrium patterns using the R code in Online Appendix B. To check 
the effects of genetic drift on these predictions, I also modeled a system of finite populations. At each generation in each deme, I 
sampled N adults from the focal deme with probability 1 — m and from each neighboring deme with probability m/2. Then I sampled 
N offspring from the frequency distribution described by Equation (9) given the adult genotype frequencies. 
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Epistasis and gene-culture coevolution: cattle as symbionts 

I implemented the Feldman and Cavalli-Sforza (FCS) model as follows (Feldman and Cavalli-Sforza, 1989). The first step is mating 
and mother-offspring transmission of the symbiote with probability v (L7 vs. u represent presence or absence of the symbiote), such 
that the offspring gene-culture distribution is 

^CUAA = ^PA ^^CUAA + -^CUAa^ 

^CUAa = ^(PaXcU^^ + pflXcUAA + -^CUAa) 
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^cuaa = pa ^^cuaa + -^cuAa^ + (1 — v) pa ^XcUaa + -^cUAa^ 

where pA and pa are the host allele frequencies. In the second step, following Feldman and Cavalli-Sforza (1989), offspring have a 
chance of picking up the symbiote via horizontal transmission. That is, each individual with trait u has probability ^[/ of being 
converted to trait U (pu is the frequency of the symbiote U among other hosts in the population). The third step is selection 
according to 



(10) 



W.UAA = 1 + Si W.uAA = 1 

W.UAa=l+Sl W.uAa=l (11) 

W.Uaa = 1 - S2 W^.uaa = 1 

R code implementing the FCS model is provided in Online Appendix C. 
Hitchhiking effects 

R code for illustrating hitchhiking effects is provided in Online Appendices D and E. These give examples of commands for iterating 
the FCS model in a series of stepping stones after secondary contact and keeping track of allele and symbiont frequencies. 
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